Thierry Phénix & Ladislas Nalborczyk
11 juin 2019
Cours 1: Introduction à l'inférence bayésienne
Cours 2: Modèle beta-binomial
Cours 3 : Modèle linéaire à 1 prédicteur continu
Cours 4 : Modèle linéaire multivarié
Cours 5: Markov Chain Monte Carlo
Cours 6: Modèle linéaire généralisé
Cours 7: Comparaison de modèles
Cours 8: Prédicteurs catégoriels et Interactions
Cours 9: Modèles multi-niveaux
Cours 10: Data Hackaton
On s'intéresse aux différences de taille entre hommes et femmes. On va mesurer 100 femmes et 100 hommes.
N = 100
men <- rnorm(N, 175, 10) # 100 tailles d'hommes
women <- rnorm(N, 170, 10) # 100 tailles de femmes
t.test(men, women)
Welch Two Sample t-test
data: men and women
t = 2.1475, df = 197.71, p-value = 0.03297
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
0.2655634 6.2364837
sample estimates:
mean of x mean of y
173.0157 169.7646
tibble(
value = c(men, women),
gender = c(rep("men", N), rep("women", N))
) %>%
ggplot(aes(x = value, fill = gender) ) +
geom_histogram() +
theme_bw(base_size = 20)
On va simuler des t-valeurs issues de données générées selon l'hypothèse de non-différence entre hommes et femmes.
nSims <- 1e4 # number of simulations
t <- rep(NA, nSims) # initialising an empty vector
t <- replicate(
nSims,
t.test(
rnorm(N, 170, 10),
rnorm(N, 170, 10)
)$statistic
)
La valeur d'un T-test est donnée par la formule :
\[ t = \frac{\mu_W - \mu_M}{S} \]
Où \( S \) est une fonction de la taille des échantillons et des variances des données.
data.frame(t = t) %>%
ggplot(aes(x = t) ) +
geom_histogram() +
theme_bw(base_size = 25)
Density, distribution function, quantile function and random generation for the t distribution with df degrees of freedom (and optional non-centrality parameter ncp).
Usage
dt(x, df, ncp, log = FALSE)
pt(q, df, ncp, lower.tail = TRUE, log.p = FALSE)
qt(p, df, ncp, lower.tail = TRUE, log.p = FALSE)
rt(n, df, ncp)
df_obs = t.test(men, women)$parameter
data.frame(t = t) %>%
ggplot(aes(x = t)) +
geom_histogram(aes(y = stat(density)),
alpha = 0.4) +
stat_function(fun = dt,
args = list(df = df_obs), size = 1.5) +
theme_bw(base_size = 25)
alpha <- .05
abs(qt(alpha / 2, df = df_obs) ) # two-sided critical t-value
[1] 1.972035
tobs <- t.test(men, women)$statistic # observed t-value
tobs %>% as.numeric
[1] 2.147452
A p-value (in black) gives the probability of observing the data we observed (or more extreme data), given that the null hypothesis is true.
t.test(men, women)$p.value
[1] 0.03297456
tvalue <- abs(t.test(men, women)$statistic)
2 * integrate(dt, tvalue, Inf, df = df_obs)$value
[1] 0.03297456
2 * (1 - pt(abs(t.test(men, women)$statistic),df_obs) )
t
0.03297456
Cours 1: Introduction à l'inférence bayésienne
Cours 2: Modèle beta-binomial
Cours 3 : Modèle linéaire à 1 prédicteur continu
Cours 4 : Modèle linéaire multivarié
Cours 5: Markov Chain Monte Carlo
Cours 6: Modèle linéaire généralisé
Cours 7: Comparaison de modèles
Cours 8: Prédicteurs catégoriels et Interactions
Cours 9: Modèles multi-niveaux
Cours 10: Data Hackaton
Comparaison de deux modèles:
\[ \underbrace{\dfrac{p(H_{0}|D)}{p(H_{1}|D)}}_{posterior\ odds} = \underbrace{\dfrac{p(D|H_{0})}{p(D|H_{1})}}_{Bayes\ factor} \times \underbrace{\dfrac{p(H_{0})}{p(H_{1})}}_{prior\ odds} \]
\[ \text{evidence}\ = p(D|H) = \int p(\theta|H) p(D|\theta,H) \text{d}\theta \]
L'évidence en faveur d'un modèle correspond à la probabilité marginale d'un modèle (le dénominateur du théorème de Bayes), c'est à dire à la likelihood moyennée sur le prior… Ce qui fait du Bayes Factor un ratio de likelihood pondéré (par le prior).
On lance une pièce 100 fois et on essaye d'estimer la probabilité \( \theta \) (le biais de la pièce) d'obtenir face. On compare deux modèles qui diffèrent par leur prior sur \( \theta \).
\[ \begin{align} \mathcal{M_{1}} : y_{i} &\sim \mathrm{Binomial}(n, \theta) \\ \theta &\sim \mathrm{Beta}(6, 10) \\ \end{align} \]
\[ \begin{align} \mathcal{M_{2}} : y_{i} &\sim \mathrm{Binomial}(n, \theta) \\ \theta &\sim \mathrm{Beta}(20, 12) \\ \end{align} \]
\[ \begin{align} \mathcal{M_{1}} : y_{i} &\sim \mathrm{Binomial}(n, \theta) \\ \theta &\sim \mathrm{Beta}(6, 10) \\ \end{align} \]
\[ \begin{align} \mathcal{M_{2}} : y_{i} &\sim \mathrm{Binomial}(n, \theta) \\ \theta &\sim \mathrm{Beta}(20, 12) \\ \end{align} \]
\[ \text{BF}_{12} = \dfrac{p(D|H_{1})}{p(D|H_{2})} = \dfrac{\int p(\theta|H_{1}) p(D|\theta,H_{1}) \text{d}\theta}{\int p(\theta|H_{2}) p(D|\theta,H_{2}) \text{d}\theta} = \dfrac{\int \mathrm{Binomial}(n, \theta)\mathrm{Beta}(6, 10)\text{d}\theta}{\int \mathrm{Binomial}(n, \theta)\mathrm{Beta}(20, 12) \text{d}\theta} \]